function [t,ux0,k,omega,A] = gen_std_wave(amp0,vr,a,b,Td)
%genRayleighStdWave Generate standard input wave for Rayleigh problem
%   Detailed explanation goes here
if nargin==4
    Td = 4.0;
end
% -----------------------------------------
%  generate wave for u0
%
fs = 200;       % sampling frequency
dt = 1.0/fs;    % delta t
% N = 800;     % number of sampling point
N = Td/dt;
t = (1:N)*dt;   % time series
t = t';
% Td = (N-1)*dt;  % duration

gaussmf = @(x, sigma, c) exp(-(x-c).^2/(2*(sigma^2)));

envelope0 = gaussmf(t,0.15,0.5);
% amp0  = 0.1;
freq0 = 3.0;
omega0 = 2.0*pi*freq0;
ux0 = envelope0.*cos(omega0*t);  % horizontal motion at (0,0)
[uxmax, maxpos] = absMaxPos(ux0);
ux0 = amp0*ux0/ux0(maxpos);

% plot(t,ux0); grid on;

domega = pi/((N+1)*dt);	% delta omega
omega = (1:N)*domega;   % omega series
omega = omega';
k = omega/vr;           % wave number

%
coef = idst(ux0);   % discrete sine transform 
%
% [coef_max,coef_max_pos] = absMaxPos(coef);
% omegad = omega(coef_max_pos);
% lr = vr*(2*pi)/omegad;  % rayleigh wave length (dominate)

%
z = 0;
eakz = exp(-a*k*z);
ebkz = exp(-b*k*z);
f1 = k.*(eakz-2*a*b/(1+b^2)*ebkz);
A = coef./f1;

end



